Nonlinear wave-wave interactions in stratified flows: Direct numerical simulations 
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To investigate the formation mechanism of energy spectra of internal waves in the oceans, direct 
numerical simulations are performed. The simulations are based on the reduced dynamical equations 
of rotating stratified turbulence. In the reduced dynamical equations only wave modes are retained, 
and vortices and horizontally uniform vertical shears are excluded. Despite the simplifications, our 
simulations reproduce some key features of oceanic internal- wave spectra: accumulation of energy at 
near-inertial waves and realistic frequency and horizontal wavenumber dependencies. Furthermore, 
we provide evidence that formation of the energy spectra in the inertial subrange is dominated by 
scale-separated interactions with the near-inertial waves. These findings support observationally- 
based intuition that spectral energy density of internal waves is the result of predominantly wave- 
wave interactions. 

PACS numbers: 47.35.Bb 
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I. INTRODUCTION 

Oceanic internal waves are the waves whose restoring 
force is buoyancy in stratified fluid. These waves are 
excited by flows over topography, tides and atmospheric 
disturbances. The energy of the waves is then transferred 
by nonlinear interactions among wavenumbers from large 
scales to small scales, and is dissipated in the small spa- 
tial scales by wave breaking. Internal waves play a signif- 
icant role in the general circulation of oceans and hence 
the climate of the Earth. 

Energy spectra of internal waves are vast with horizon- 
tal wavelengths varying from 10 meters to 10'°^ meters, 
vertical wavelengths from 10 meters to 10'^ meters, and 
time periods from 10^ seconds to 10^ seconds. 

The complexity of the internal-wave fields arises not 
only from its extended range of scales, but also from their 
interactions with the other major players in ocean dy- 
namics including eddies, mean currents and shear flows. 
The main dynamical role of the internal waves is to store 
energy and transfer it across different scales and large dis- 
tances. Hence the waves constitute a large and complex 
geosystem containing a broad range of interacting scales 
and affecting significantly most of the active players in 
ocean dynamics. 

However, surprisingly and despite all the complexi- 
ties, energy spectra of the internal waves in the oceans 
appear to be somewhat universal. It is given by the 
Garrett-Munk (GM) spectruniiiSJ,. It is believed that 
the internal-wave spectrum may be a result predomi- 
nantly, if not exclusively, of nonlinear interactions among 
waves. 

In this paper we test the hypothesis by direct numer- 
ical simulations of the reduced model of stratified rotat- 
ing turbulence. Our model contains only wave modes 
of stratified wave turbulence, and completely excludes 
vortices and horizontally uniform vertical shears. De- 



spite these simplifications, our numerical model repro- 
duces some key features of spectral energy density ob- 
served in the oceans. 

As a result of appearance of the GM spectrum, internal 
waves in the ocean have been a subject of active research 
ever since. An important milestone was a review by 
Miiller et al.— . It focuses on the resonant wave- wave in- 
teraction theory, called the weak turbulence theory. The 
result of the review is that the resonant wave-wave inter- 
actions in the stratified wave turbulence are dominated 
by specific, "named" nonlocal wave-wave interactions in 
the wavenumber space. The classification of the "named" 
nonlocal interactions appeared first in McComas-. 

Since then it was understood that in the oceans local 
spectral energy density may change owing to the propa- 
gation of energy from other parts of the ocean as well 
as owing to the nonlinear wave-wave interactions. A 
useful and intuitive diagram in the wavenumber space 
that separates the two regions appear in^. More re- 
cently, Levinei proposed modification of the GM spectra 
to take into account the dependence of the characteristic 
depth as a function of frequency lu. Furthermore, Lvov 
and Tabak— developed a novel Hamiltonian structure for 
waves in stratified flows that we use below for our nu- 
merical modeling. Historical observational data were re- 
viewed in Lvov et al.— , where major deviations from the 
GM spectrum were categorized. Useful phenomenologi- 
cal characterization of fluxes of energy in internal waves 
were put forward by Polzin—. 

However, the GM spectrum stood the test of time and 
still stands as the canonical model of the internal-wave 
spectral energy density. The GM spectrum is a func- 
tion of the frequency uj and the vertical wavenumber m, 
E{uj,m). However, it is rather hard and expensive to 
measure the spectrum as a function of both uj and m. 
At least in the 1970's, when the series of the GM spec- 
tra are publishedi'^'^, only one-dimensional spectra were 
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generally available (with the exception of the IWEX ex- 
periment). In particular, i?timc(w) was obtained from 
time series of mooring current meters and i?vcrticai('7i) 
from vertical profilers. Furthermore the large wavenum- 
ber limit of -Bvorticai(w) has dependency. 

Garrett and Munk proposed that the spectrum is sepa- 
rable, that is the product of a function of uj and a function 
of m: 

E{UJ, m) OC i^timo('^)£^vortical("l)- (1) 

Then £'timo(w) and -Everticai(w) were properly normal- 
ized and chosen so that the resulting spectrum matches 
experimentally measured one-dimensional oj and m spec- 
tra. The resulting spectrum is given by Eq. (|22p below. 

The assumption of the separability allows one to con- 
struct a function of two arguments, E[uj, m) out of two 
one-dimensional functions, Etimd'-^) and vertical (™)- 
However, if one relaxes the assumption, there is more 
than one way to obtain the two-dimensional spectrum, 
E{uj,m), to fit the observed one-dimensional spectra, 
-Btime(w) and -Everticai(TO)- Morcover, it recently became 
quite apparent that the assumption ([T]) is not satisfied 
in the oceansii. Our numerical simulations also demon- 
strate that the assumption U]) is not satisfied uniformly. 
It also appears in our direct numerical simulations that 
resulting energy spectra do not display universal behav- 
ior. Rather, our simulations demonstrate accumulation 
of energy around the horizontally longest waves. Further- 
more, we argue based on our direct numerical simulations 
that the energy spectrum in the inertial subrange is de- 
termined by nonlocal interactions with the accumulation. 
The nonlocal interaction in the wavenumber space is one 
of the "named" nonlocal interactions identified by Mc- 
Comas^. As a result of the nonlocal interactions, the 
details of behavior of the stratified wave turbulent sys- 
tem depends upon details of the accumulation. Conse- 
quently the resulting energy spectra are non-universal in 
our numerical experiments. 

Despite apparent non-universality and non- 
separability of spectral energy density in the inertial 
subrange, our simulations do exhibit certain key features 
that are observed in the ocean. Namely, our simulations 
have clear peaks at inertial frequencies which correspond 
to the accumulation of energy, as observed in the ocean. 
Our largest numerical simulation does demonstrate the 
cj"^ dependence of the energy spectrum that appears 
prominently in moored observations. Furthermore our 
simulation demonstrate realistic fc~^ dependence on the 
horizontal wavenumbers. The behavior of the spectra 
in the inertial subrange, formation of accumulation, 
apparent non-universality and violation of separability 
([T|) can be qualitatively interpreted in terms of the 
nonlocal "named" interactions. These findings support 
observationally-based intuition that spectral energy 
density of internal waves is formed primarily by the 
nonlinear wave- wave interactions. 

The stratified rotating wave turbulence is a subset of 
a much more complicated subject of rotating stratified 



strong turbulence governed by the Navier-Stokes equa- 
tion. Complexities of the rotating stratified strong tur- 
bulence appear owing to coexistence of waves, shears 
and vortices and their interactions. The rotating and 
stratified strong turbulence has been a subject of inten- 
sive research in last few decades. Extensive numerical 
studies of rotating turbulence, stratified turbulence and 
turbulence with both rotations and stratification were 
performedi^iii^ii^ii^. Accumulation of energy at the hor- 
izontally largest scales were reported in direct numerical 
simulationsi^ii^. The accumulation in the rotating strong 
turbulence happens presumably owing to the inverse cas- 
cade of two-dimensional turbulence. In contrast, in the 
stratified rotating turbulence, the resonant wave-wave in- 
teractions commonly cause the accumulation of energy at 
the horizontal largest scales. 

We emphasize that the rotating stratified wave turbu- 
lence is dominated by nonlocal wave-wave interactions 
in the wavenumber space. Anisotropic wave turbulence 
systems often exhibit nonlocal interactions: drift wave 
turbulence, Rossby waves and MHD turbulenc o^^'^^ . 
The scenario is in contrast to the local interactions of 
isotropic wave turbulence systems. In the locally inter- 
acting systems the spectra are insensitive to the details 
of large-scale and small-scale motions, and are the re- 
sults of the local interactions among wavenumbers. Ex- 
amples of the locally interacting systems include waves 
on water surfaces. In particular, universal behavior 
was observed in direct numerical simulations of gravity- 
wave and capillary-wave systemsi?'— i2°'lii??.. Note that 
isotropic Navier-Stokes turbulence is also widely believed 
to be a locally interacting system. 

The paper is written as follows. In Sec. [ITl we give 
the detailed description of our numerical model and as- 
sumptions used. In Sec. IIIII we elaborate on our numeri- 
cal methods, and explain pumping and damping mecha- 
nisms. In Sec. llVl we account for results of our numerical 
experiments. The formation mechanism of energy spec- 
tra is discussed in Sec.|Vl Section IVll provides summary. 



II. HAMILTONIAN FORMALISM FOR 
INTERNAL WAVES 

In this section we provide a description of the model 
that we use for our numerical study. The model is based 
on the Hamiltonian description of the wave modes of 
the incompressible stratified rotating flows in hydrostatic 
balance. The Hamiltonian description appeared in Lvov 
and Tabak-§ and is presented here for completeness. Our 
model explicitly excludes vortices and horizontally uni- 
form shear flows. Despite the simplification, the result- 
ing spectrum does display some key features that are ob- 
served in the ocean. 

As a starting point, we take the equations of motion 
satisfied by an incompressible stratified rotating flow in 
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hydrostatic balance under the Boussinesq approximation: 



d dz 
di'd'p 
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(2) 



These equations are derived from mass and horizontal- 
momentum conservations and hydrostatic balance. The 
equations are written in isopycnal coordinates with the 
density p replacing the height z in its role as indepen- 
dent vertical variable. Here u = {u, v) is the horizontal 
component of the velocity field, and — {—v,u). The 
gradient operator V — {d/dx, d/dy) acts along isopyc- 
nals. The inertial frequency due to the rotation of the 
Earth / is assumed to be constant, g is the acceleration 
of gravity, and po is a reference density (in its role as in- 
ertia) which is taken to be a constant in the Boussinesq 
approximation. Finally M is the Montgomery potential 



M = P 



gpz. 



The expression for the potential vorticity in these co- 
ordinates is, 



9 = 



/ + v- 
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(3) 



where = {—d/dy, d/dx) is the two-dimensional rota- 
tion operator. Here we introduced 

j^_pd^ _ dz_ 
g dp"^ ^ dp 

to be a normalized differential layer thickness. The po- 
tential vorticity is conserved along particle trajectories. 
Since the fluid density is also conserved along the trajec- 
tories, an initial profile where the potential vorticity is 
a function of the density will be preserved by the flow. 
Hence any initial chosen profile will stay in the fluid all 
the time. This observation allows us to propose that 
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Here lioip) is a reference stratification profile defined as 



no(p) = - 



N{py 



(5) 



and N{p) is the buoyancy (Brunt -Vaisala) frequency, 
which we shall regard as a constant Nq. 

In order to separate wave and vorticity dynamics, we 
decompose the fluid velocity into its gradient and rota- 
tional parts, i.e. 



(6) 



The vorticity is derived from the potential vorticity and 
is distinct from the vorticity in Cartesian coordinates. In 
terms of the potentials </> and -0, the constrain (j4]) reads 

/ -f AtA = goH . 

Therefore we can express ill as a. function of H so that 
Eqs. @ and ^ are satisfied. As a result, if wc redefine 
n as n + Ho, Eqs. ([2]) reduce to the pair: 



dW 

— +V - ((H + Ho)(V(/. 



V^A-i(gon-/))) =0, 



v^A-i(gon-/)| 



-A^V 
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n(pi) 
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vA-i(gon-/))) 



dpidp2 = . 



(7) 



The pair of equations ([7]) form a canonical conjugate 
pair of the Hamiltonian equations, 

dn _ sn dcj) _ sn 

'dt ~ ^~5^ ' 'dt ~m' 

The Hamiltonian is the sum of kinetic and potential en- 
ergies, 
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(8) 



n = I dxdp I 



p' 



{n{x,p)+no)\\/cb{x,p)+W^A-\oU{x,p)\ .(9) 



We therefore consider in this paper the following re- 
duced model with the following assumptions and con- 
straints: 

• only wave motions are considered 

• vortex motions are excluded 

• horizontally uniform vertical shear is excluded 

• potential vorticity is constrained to be constant on 
an isopycnal surface. Its value is determined by the 
underlying rotations 

• constant underlying rotation / = const is assumed. 
Thus we explicitly exclude the /3 effect. 

• hydrostatic balance is assumed. 

• buoyancy frequency N{z) = const is assumed to be 
constant with depth. 

We do realize that the oceans are more complicated 
than these idealizations. Certainly for general ocean 
modeling the above assumptions constitutes a gross sim- 
plification. The reason for our choice is that we would 
like to find out whether it is sufficient to study wave- wave 
interactions alone to determine the form of the internal- 
wave spectral energy density. We show below that this 
is indeed the case to a large extent. We show that our 
reduced model reproduces key characteristic behavior of 
the oceans. The effects of relaxing the above assumptions 
is the subject of future work. 
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III. DETAILS OF NUMERICAL METHODS: 
IMPLEMENTATION AND INTERPRETATION 

A. Wave turbulence 

To proceed, we perform Fourier transformation and 
canonical transformation to the field variable, a{p), de- 
fined as 



a{p) 




n(p) 



9_\k[-^ 
2uj No' 



[P) 



(10) 



with linear coupling of the Fourier components of the 
stratification profile, 11, and the horizontal velocity po- 
tential, (/>. The three-dimensional wavenumber, p, con- 
sists of a two-dimensional horizontal wavenumber in the 
isopycnal surface, k, and a vertical density wavenumber, 
m. The linear frequency in the isopycnal coordinates is 
given by the dispersion relation. 



9" |fcp 



(11) 



The usual vertical wavenumber, /c^, and the density 
wavenumber, m, are related as to = —g / [pQNQ)kz. 

Then, the pair of canonical equations of motion ([8| is 
rewritten as a single canonical equation. 



.da{p) Sn 



dt 5a* {p) 



(12) 



with the standard Hamiltonian of three-wave interacting 
systems--^ , 

j dpLo{p)\a{p)\^ 

dpdpidp2 {{V^^^p^a{p)a*{pi)a*{p2) + c.c.) 

+ iUp,pi,p2a{p)aipi)a{p2) + c.c.)) • 

(13) 

Here, 5 /5a* is the functional derivative with respect to 
a*{p), which is the complex conjugate of a{p), and the 
abbreviation c.c. denotes complex conjugates. The ma- 
trix elements, V^^ and Up_p^^p^, have exchange symme- 



tries such that V^^ p^ 



^P2,Pl ^P,Pl,P2 



P1,P,P2 ■ 

The Hamiltonian is the canonical form of the 

Hamiltonian of wave turbulence system dominated by 
three-wave interactions^^. The first term describes lin- 
ear noninteracting waves, the second term correspond to 
the nonlinear three-wave scattering processes. The wave 
turbulence theory provides a powerful framework to de- 
scribe spectral energy transfer in the systems dominated 
by wave- wave interactions. A detailed review of the wave 
turbulence theory and its applications to internal waves 
is outside of the scope of the present paper, and is given 
in Lvov et al.— . Here it is sufficient to note that there 
are the following important classes of scale-separated res- 
onant interactions among waves^: 



• The vertical backscattering of a high frequency 
wave by a low frequency wave of twice the vertical 
wavenumber into a second high frequency wave of 
oppositely signed vertical wavenumber. This type 
of scattering is called Elastic Scattering (ES). 

• The scattering of a high frequency wave by a low 
frequency and small vertical wavenumber wave into 
a second, nearly identical high frequency and large 
vertical wavenumber wave. This type of scattering 
is called Induced Diffusion (ID). 

• The decay of a small wavenumber wave into two 
large vertical wavenumber waves of approximately 
one-half the frequency. This is called Parametric 
Subharmonic Instability (PSI). 

The classification provides a useful interpretive frame- 
work to characterize resonant wave-wave interactions in 
stratified flows. We will show below that results of our 
numerical simulations can be qualitatively interpreted by 
using this classification. Detailed theoretical analysis of 
scale-separated interactions of this type will be presented 
in Lvov et al.— . 



B. Numerical setting 

To achieve non-equilibrium statistically (near-)steady 
states, we have to model both processes of pumping en- 
ergy to the internal-wave field and of damping energy 
from the field. The processes of the pumping include in- 
teractions with surface waves and tides. How to model 
these processes in the wavenumber space is the subject 
of present oceanographic research. We model the pump- 
ing processes phenomenologically. In what follows we 
assume that the pumping occurs on large length scales, 
and is relatively local in the wavenumber space. The 
processes that remove energy from the wave field include 
wave breaking, turbulent dissipation, reflection from sur- 
face and bottom boundary layers and interaction with 
topography. Again, the spectral details of the processes 
is a subject of intensive research. We assume that the 
processes are especially effective for small length scales 
(large wavenumbers). 

In our numerical simulations, most of the energy pro- 
vided by the external forcing accumulates around the 
horizontally longest waves owing to nonlinear interac- 
tions among waves. The horizontally longest waves have 
frequencies near the inertial frequency /. Therefore the 
waves are called "near-inertial waves." Then the energy 
is transferred through the inertial subrange, where there 
is no signiflcant pumping or damping. Note that the in- 
ertial subrange therefore refers to the range in wavenum- 
ber space without effective forcing and damping, while 
the near-inertial waves refers to the waves with the near- 
inertial frequency. Subsequently energy is absorbed in 
the dissipation range. It is the nonlinear interactions 
among internal waves that determine the form of the 
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TABLE I: Numerical parameters. The wavenumbers, k and m, are discretized and they have integer values. 



modes / (xlO ''rad/sec) Lv (xlOkg/m^) 



V2/3 
0.25 
1 

\/2/3 
72/3 



forcing 



initial condition 



I 

II 

III 
IV 



512^ X 256 

256^ 

256^ 
512^ X 256 



V 102r X 512 



2.7 
5 
5 

2.7 
2.7 



|fcp 



none 

+ < 6^ 



none 
a; -3/ 



GM 
white noise 
white noise 
GM without long waves 
white noise 



spectrum in the inertial subrange and the formation of 
accumulation of energy at the near-incrtial waves. The 
nonlinear interactions among waves is the main focus of 
the present paper. 

In the direct numerical simulations, we add external 
forcing and hyper-viscosity to the canonical equation 
(|12p . Thus the dynamic equation used in the simulations 
is given by 

^ - -iuj{p)a{p)+N{a{p)) + F{p) - D{p)a{p){U) 

Details of numerical algorithm are the following: The 
linear terms, i.e. a linear dispersion term and a dissi- 
pation term, ~D{p)a{p), are explicitly calculated. The 
nonlinear terms, Af{a{p)), are derived from the nonlin- 
ear parts of the canonical equation p2p with Hamiltonian 
(lisp . They are obtained numerically by a pseudo-spectral 
method with the phase shift under the periodic boundary 
conditions for all three directions. The external forcing, 
F{p), is implemented by fixing the amplitudes of several 
small wavenumbers to be constant in time. This is com- 
monly used to simulate forcing in numerical experiments. 
The dissipation is modeled as hyper- viscosity: 

D{p) = Di,\k\'' + D^\m\\ (15) 

Here, Z?h and are chosen so that the dissipation is ef- 
fective for wavenumbers larger than the half of maximum 
wavenumbers. 

The wavenumbers are discretized as p = 
(27r/Lhfc, 27r/Lvm), where Lh and Lv are horizon- 
tal periodic length and vertical period in the isopycnal 
coordinates, and k and m are integer-valued wavenum- 
bers. We are going to use integer-valued (dimensionless) 
wavenumbers from here on in this paper. Time-stepping 
is implemented with the fourth-order Runge-Kutta 
method. In all the simulations, the buoyancy frequency 
and horizontal period are fixed at A'o = lO^^rad/sec 
and = lO^m, respectively. 

We perform a series of five numerical experiments that 
are listed in Table [J The total energy per unit peri- 
odic box of all the numerical experiments except Run 
V is around 3 x 10'^ J/ (kg • m^) which is characteristic 
of the oceans. The total energy in Run V is around 
1.2 X 10'^ J/ (kg • m^). The values of total energy density 



and the dissipation rate of Run V in Cartesian coordinate 
are 1.3 x lO^'^J/kg and 5.0 x 10~^^W/kg, respectively. 
These values are in good agreement with observations 
and theoriesi^. 

Stratified rotating turbulence can be characterized by 
two dimensionless numbers, the Rossby number Ro, 
which is the ratio of the inertia force to the Coriolis 
force, and Richardson number Ri, which is the ratio of 
the buoyancy to the shear. They are defined as 

Ro - {\{u-\/)u\/\fu^\), 
Ri = {N^/S^), 

where (•) denotes averaging in the numerical box and N 
and S show local buoyancy frequency given by local strat- 
ification and local shear, respectively. The values of the 
Rossby number and the Richardson number measured in 
Run V are 2.7 x 10^^ and 8.6 x 10^, respectively. 

C. Integrated and cross-sectional spectral energy 
density 

The present paper is concentrated on the investiga- 
tion of the behavior of spectral energy density in the 
inertial subrange. As mentioned above our numerical 
scheme employs the pseudo-spectral algorithm. Con- 
sequently it is convenient to concentrate our attention 
on the (fc,m) spectra. Indeed, the numerical box is a 
triple periodic box in the wavenumber space. On the 
other hand, oceanographers find it convenient to think 
in (cj, m) space. Indeed, the to and m spectra could be 
measured experimentally. 

There are merits to both ways of thinking. For exam- 
ple, we will demonstrate below that the spectrum is more 
separable in (fc, m) space for Run II and III. On the other 
hand. Run V could be better interpreted in {lu, m) space. 
The oceanic spectra also tend to be more separable in 
(tj,m) space. Indeed, the GM spectrum is assumed to 
be separable in (cj, m) space. We will use both ways of 
thinking through the paper. 

Two-dimensional energy spectra are measured in the 
shell of radius k as 

E{k,\m\)^ ^ ^ c.|a(fc',sm)|2. (16) 

k-l/2<\k'\<k+l/2 s=±l 
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Integrated energy spectra are defined as 



i?i„t(fc)=5]i?(fc,|m|), (17) 

m 



and 



^int(|m|) = ^i?(fc,|™|). 



(18) 



Cross-sectional spectra, Em{k) andi?fc(|m|), are obtained 
from the two-dimensional energy spectra as a function 
of horizontal wavenumbers k along a certain density 
wavenumber m and as a function of density wavenumbers 
m along a certain horizontal wavenumber k, respectively. 

As well as in (k,m) space, we numerically obtain en- 
ergy spectra in (w, m) space. It is made in the similar 
way how observational spectra are. Namely, we choose a 
point on the "surface" of our numerical ocean. We then 
conduct a "vertical mooring," recording a time series of 
the horizontal velocity u{xQ,p;t) at a fixed horizontal 
position Xq. The kinetic energy spectra is defined as 



K{uj,m) = ^\u{xo,m;uj)\'^, 



(19) 



where u{xQ,m]uj) is the Fourier component of the hor- 
izontal velocity with respect to the vertical and time 
series. The integrated kinetic energy spectra K-mti^) 
and ifintdTij), and the cross-sectional kinetic energy 
spectra Km{uj) and ifi^dml) are defined from the two- 
dimensional kinetic energy spectrum, K {10,111,), in the 
same way as the energy spectra in (k,m) space. 

As explained in the introduction, in the 1970's only 
one-dimensional spectra were widely available with the 
exception of the IWEX experiment. To measure 
EvcTticaii'm) spectrum in the ocean, one can use vertical 
profilers at a given position. To measure i?timc(^) spec- 
trum, one calculates time series of the mooring current 
meters. Garrett and Munk assumed that the spectrum is 
separable in (uj,m) space. In other words, the spectrum 
is the product of function of lu alone and function of m 
alone. In functional form, this statement is written as 
Eq. _ 

Then i?timo('-^) and i?vorticai(™) were properly normal- 
ized to reproduce the characteristic total energy density 
of internal waves. Then these functions were chosen so 
that the resulting (w, m) spectrum is consistent with both 
LU spectrum and m spectrum. In particular, moored spec- 
trum was chosen to be 



1 



(20) 



Notice that it has an integrable peak at the inertial fre- 
quency /. Moreover it has w"^ dependence for high fre- 
quencies, as is prominently displayed in moored observa- 
tions. Then the m spectrum was chosen as 



irtical 



(m) OC 



1 



(21) 



Here m* is the characteristic wavenumber determined by 
scale height. Detailed analysis of these choices will be 
presented in Polzin and Lvov—. 

The assumption of separability ([1]) allows one to con- 
struct a function of two arguments, E{ll!, m) out of two 
one-dimensional functions. The resulting GM spectrum 
in {uj, ni) space is given by 



1 



1 



(22) 



However, if one relaxes the assumption of the separa- 
bility ID), it is more than one way to obtain the two- 
dimensional spectrum, E{uj,m), to fit the observed one- 
dimensional spectra, Eqs. (PO)) and (HI]). Moreover, it 
recently became quite apparent that the separability ([1]) 
is not satisfied in the oceans^i. Therefore, we could ar- 
gue that reality and accuracy of the separability ([1]) per- 
haps deserves a closer examination. Furthermore, it may 
be possible that the one-dimensional spectra admits ex- 
planation without the separability. We therefore stress 
through the paper the conceptual difference between the 
integrated Eiat{uj) and Eint{m) spectra that could be 
confirmed by oceanographic measurements and the two- 
dimensional spectra which is a function of both lu and 

TO. 

A detailed analysis of the observational data is outside 
of the scope of the present paper. We will present such 
analysis in Polzin and LvovSl. There details of IWEX 
and other experiments will be presented and reanalyzed. 
We will also present there results and analysis of more 
modern observations. Brief catalog of historically avail- 
able observational programs in the past three decades, 
along with characterization of deviations from the GM 
spectrum, is available in Lvov et al.^. 

The linear dispersion relation (jlip links lu and m 
with the horizontal wavenumber k. With the linear 
dispersion relation, the spectrum can be transformed 
from both wavenumber space, [k,m), into frequency- 
horizontal- wavenumber space, [k,Lu), or frequency- 
vertical- wavenumber space, (m,w). In particular, the 
GM spectrum in (k, m) space is given as 

du 



EQM{k,m) = £;gm(w, to) — 



1 



1 



OC 



2 J 9!_Mi ItoI (to^ + TO*^ 



f 



- (23) 



By using Eq. (|23p we can define the integrated to spectra. 



E 



int,GM 



(to) 



No 



f 



EGM{uJ,'rn)dLU 



I 



EGM{k,m)dk.{2A) 



In the large- wavenumber and large- frequency limit, 
Eq. has the self-similar asymptotic form given by 



E GM {^^ , rn) OC LU ^|to|" 



(25) 
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Similarly, in the large- wavenumber limit, Eq. (|23p has a 
form 

EGM{k,m)(xk-^\m\~\ (26) 

Both Eqs. (|25p and (|26|) should give the same one- 
dimensional spectrum as a function of vertical wavenum- 
bers, 

£^int,GM(™) « |m|~^. (27) 

The difference between the exponents of m in Eqs. (|26|) 
and (P7)) comes from the non-separability of Eq. (|23p in 
the (fc, to) space. 

Note that the power-law exponents of one-dimensional 
spectra and those of two-dimensional spectra do not al- 
ways coincide. This could be explained by the fact that 
the exponents of the integrated spectra are a weighted 
mean of those of the cross-sectional spectra, see (fT7|) and 
(fTS]) . Therefore, the exponents of the integrated spectra 
can be different from those of the cross-sectional spectra 
unless the spectra are separable. 

IV. NUMERICAL RESULTS 
A. Run I 

Since the GM spectrum was believed to be the uni- 
versal oceanographic spectrum, the natural question to 
ask is whether the GM spectrum could be a steady state 
of our numerical wave model. We realize that such a 
guess may be far-fetched, but it is nevertheless appeal- 
ing to check it. Therefore we are led to examine sta- 
tistical stability of the GM spectrum. We assume that 
the steady state of the ocean can be characterized by 
an energy flux from pumping regions to damping regions 
in the wavenumber space. Therefore to achieve a truly 
steady state we need to model both pumping processes 
and damping processes. As explained above, in modeling 
the damping we use traditional hyper-viscosity approach 
in Eq. (fT5|) . As we will see below, the resulting spec- 
trum depends upon the large-scale flows. Therefore, the 
only way to see whether the GM spectrum is close to 
being a steady state that is independent of the form of 
the pumping is to choose no pumping at all. Thus we 
are led towards modeling the system as a freely decaying 
system, i.e. without external forcing. 

The GM spectrum with the cut-off in large horizon- 
tal and density wavenumbers is employed as the initial 
energy spectrum. The cut-off is introduced to avoid de- 
creasing accuracy of the pseudo-spectral method. The 
initial phases of complex amplitude, a(p), are given by 
uniformly distributed random numbers in [0, 27r). 

Figure [T] shows the initial spectrum and the energy 
spectrum after about 35 hours in ocean time. If the GM 
spectrum were to be a universal steady state of our nu- 
merical model, it would change very little in this simula- 
tion. Instead, the density exponents rapidly change from 
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FIG. 1: Cross-sectional spectra Ek{\m\) of GM spectrum as 
the initial condition (upper) and of energy spectrum after 
about 35 hours (bottom) of Run I. Significant difi'erences in 
the region 16 < |m| < 64 indicates that GM spectrum is not 
statistically steady. 



— 1 to —2 only after 1.5 days. This spectral change occurs 
faster than dissipation effects with timescale estimated to 
be about 50 days at |to| = 32. Therefore it appears that 
the GM spectrum is not a stable universal spectrum in 
the wavenumber region 16 < \m\ < 64 at least for this 
particular numerical experiment. The behavior of nu- 
merical experiment appears to contradict the arguments 
Qj5j26 based on induced diffusion approach that all energy 
spectra are rapidly relaxed to the GM spectrum, but is 
consistent within. 



B. Run II and III 

We certainly realize that our numerical model of Run 
I may not fully describe the ocean. Then we question 
what could be the nature of the statistical steady state 
of our wave model of stratified rotating turbulence. Here 
we model the pumping phenomenologically. We assume 
that waves are forced at horizontally and vertically large 
scales, i.e. small wavenumbers. In particular we choose 
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the forced wavenumbers such that 

Ifcpp + Impp < 6^ (28) 

The forcing is modeled by fixing the amphtude of the 
forced wavenumbers to be constant in time. We also 
added an additional external dissipation —D^{uj — Nq)'^ 
for uj > Nq in these runs to avoid violation of the hy- 
drostatic balance. Detailed investigation of the non- 
hydrostatic effects is outside of the scope of present pa- 
per. 

Since the energy flows in the wavenumber space 
strongly depend on the inertial frequencies /— , we per- 
formed Run II and III with differing values of the inertial 
frequencies / (see Table HI. For Run II, the inertial fre- 
quency / = 0.25 X 10"'* rad/sec, while for Rmi III the 
value of the inertial frequency is / = 1 x 10^"* rad/sec. 
These inertial frequencies correspond to latitudes of 10 
degrees and 45 degrees. 

Note that there is a "critical" latitude where the fre- 
quency of the semidiurnal tide is equal to twice the in- 
ertial frequency. The principal difference between these 
two runs is that Run II is southwards of "critical" latitude 
in the northern hemisphere while Run III is northwards. 

Our phenomenological forcing ([28ll corresponds to the 
band of frequencies greater than 4/ for Run II and -\/2/ 
for Run III. These values are the result of the discrete- 
ness of the numerical grid. Certainly the forcing is not 
necessarily characteristic of the ocean and is purely phe- 
nomenological, and chosen for numerical simplicity and 
to allow easy interpretation. 

The kinetic energy spectra as functions of uj and m of 
Run II and III are shown in Figs. [2] and [3l Both frequency 
spectra (left in Figs. and [3]) have peaks at the inertial 
frequencies. This indicates that most energy accumulates 
in the near-inertial frequencies. The behavior is charac- 
teristic of the ocean. Indeed, these peaks correspond to 
the integrable singularity at the inertial frequency of the 
GM spectrum 

After about 10"^ days from the initial time when all 
the wavenumbers have extremely small energy, the sys- 
tems are still transient and the accumulation of energy 
in the near-inertial frequencies has not reached complete 
statistically steady states. Thus the timescale of the de- 
velopment of the accumulation is relatively large. Con- 
sequently one may conjecture that other processes, not 
present in our numerical model may affect the ocean at 
the large timescales. The most important of such pro- 
cesses which is not included in our simulations are the 
P effects. Indeed, inclusion of (3 effects would lead to 
existence of Rossby waves in our simulation, which may 
alter our results significantly. Other effects not present 
in our simulations, which can affect larger scales, may 
also alter the results. However, our subject is not to 
see how the energy spectra in the large-scale flows are 
formed but rather to investigate how the energy spectra 
in the inertial subrange are formed. For this subject, our 
numerical model is consistent with observationally-based 
intuition that wave-wave interactions is the predominant 
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FIG. 2: Integrated kinetic energy spectra Kint{(^) and 
Ji'int(l'Til), and cross-sectional energy spectra Kmiuj) and 
J^ij(|m-|) of Run II. The kinetic energy spectra are obtained 
as functions of u and m from time and vertical series of 
the liorizontal velocity u{xo, p;t). The GM spectrum scales 
as uj~'^\m\~^ for large frequencies and density wavenumbers. 
The cross-sectional spectra that are functions of frequencies 
(left) are shown every eight curves for visibility. The cross- 
sectional spectra that are functions of vertical wavenumbers 
(right) are shown when a; = 2, 4, 8, 16, 32 x 10~'*rad/sec. 



mechanism of forming the internal- wave spectrum in the 
inertial subrange. 

The energy spectrum of Run II has only weak accumu- 
lation of energy in the near-inertial frequency. This could 
be qualitatively explained by the fact that all the forced 
wavenumbers have frequencies greater than 4/ for Run 
II. Consequently, PSI is effective in transferring energy 
to frequencies around 2/ and large vertical wavenumbers. 
Then PSI can no longer be effective in transferring energy 
towards the accumulation of energy in the near-inertial 
frequencies. Indeed, second PSI transfer would make the 
vertical wavenumber much larger and thus would move 
away from the accumulation. On the contrary, the en- 
ergy spectrum of the Run III has moderate accumulation, 
stronger than that of Run II. This could be explained by 
the fact that the forced modes have frequencies greater 
than ^/2/. Then PSI can be effective in transferring en- 



9 




1 2 4 B 1C 32 H 123 



FIG. 3: Integrated kinetic energy spectra, and cross-sectional 
energy spectra of Run III. See the caption of Fig. [2]for details. 



ergy from the wavenumbers that have frequency 2/ to 
the accumulation whose frequencies are close to /— . 

Since the forcing does not have just one frequency 
but several frequencies, the integrated spectra Kint{i^j) 
have huge oscillations corresponding to the linear fre- 
quencies of the forced wavenumbers. Therefore, the 
power-law region of the frequencies, 10~* < w < 10~^, 
are strongly contaminated by the forcing. However, the 
cross-sectional spectra with |m| > 6 are not much af- 
fected by the forcing. 

The vertical- wavenumber spectra (right in Figs. [2] and 
[3]) appear to exhibit neither separable nor self-similar 
patterns. On the other hand, the integrated spectra do 
exhibit the self-similarity. The integrated spectra of the 
vertical wavenumbers are made mainly from the peaks 
at the near-inertial frequencies. Therefore, the cross- 
sectional spectra do not always have similarities to the 
integrated spectra. 

It appears that it is advantageous to obtain {k, m) 
spectra as well. The integrated energy spectra Eint{k) 
and Eint{\m\), and cross-sectional energy spectra Em{k) 
and Ek{\m\) of Run II and Run III are shown in Fig. [4] 
and in Fig. [51 respectively. The figures indicate that the 
energy spectra as functions of k and m appear to be more 
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FIG. 4: Integrated energy spectra i?int(fc) and _Bint(|77i|), and 
cross-sectional energy spectra Emik) and i?fc(|m|) of Run II. 
The GM spectrum scales as fe~^|m|~^ for large horizontal and 
density wavenumbers. The cross-sectional spectra are shown 
every four curves for visibility. 

separable than the kinetic energy spectra as functions of 
Lo and TO in Run II and III. 

The integrated and cross-sectional spectra in Figs. 2] 
and [5] have power-law regions. The power-law exponents 
of the integrated and cross-sectional spectra are shown 
in Figs. [B] and [71 The power-law exponents of the cross- 
sectional spectra by the least-square fitting as 

£;™(A:) cx fc"(l"l) in fee [8, 40] 

and 

Ek{m) (X |to|'^('') in |to| € [8,40]. 

The power-law exponents of the integrated spectra are 
also obtained as 

Ei^t{k) oc fc"'- 

and 

;Bint(|m|) cx \mf"'^ 
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FIG. 5: Integrated and cross-sectional spectra of Run III. 
See the caption of Fig. |4]for details. 

in the same wavenumber regions. In Run II the energy 
spectrum is close to double-power 

E{k,m) (X k-^\m\-^-^ (29) 

in the large horizontal and density wavenumbers. 

The exponents of the horizontal wavenumbers of both 
integrated and cross-sectional spectra, aint and a(m), 
roughly agree with the those of the GM spectrum, —2 
of Eq. (|26)) . This fact supports that the oceanic spec- 
tra of the internal waves can be explained by the wave- 
wave interactions. The integrated spectra are the sum 
of the cross-sectional spectra. The integrated spectra of 
the horizontal wavenumbers in the inertial subrange are 
not affected by the accumulation of energy around the 
horizontally longest waves. 

Indeed, to obtain the integrated spectra as functions of 
horizontal wavenumbers, we use Eq. (|17p and integrate 
over all m value for each k value. Thus, the integrals 
over vertical wavenumbers in large horizontal wavenum- 
bers are unaffected by the strong presence of the near- 
inertial accumulation. Therefore the integrated spectra 
of k in large horizontal wavenumbers well reflect the be- 
havior of the spectra in the inertial subrange, and are 
immediately insensitive to the details of the accumula- 



FIG. 6: Power-law exponents of each cross-sectional spec- 
trum in Run II. left: a(|m|), which are the exponents of the 
cross-sectional energy spectra Em{k), right: /?(fe), which are 
the exponents of the cross-sectional energy spectra i?fc(|m|). 
The error bars are obtained by fitting errors due to the least- 
square method. The exponents of the integrated spectra are 
also shown. 



tion of energy at small horizontal wavenumbers. Similar 
arguments can be applied to the integrated spectra of fre- 
quencies. Therefore, the integrated spectra of w in high 
frequencies are determined by the behavior in the iner- 
tial subrange and are insensitive to the presence of the 
accumulation of energy at the near-inertial waves. 

Consequently, the exponents of horizontal wavenum- 
bers of both the integrated and cross-sectional spectra, 
—2, are characteristic in our numerical experiments, and 
is consistent with oceanographic observations. 

On the other hand, to obtain the integrated spectra 
as functions of vertical wavenumbers, one should use 
Eq. (|18p to integrate for each m value over all k val- 
ues. Therefore, when integrating over k values, the value 
of the integrals are always strongly affected by the ac- 
cumulation of energy at small k values. Therefore the 
exponents of the integrated spectra of m are determined 
exclusive by the near-inertial accumulation, and are in- 
sensitive to the behavior in the inertial subrange. Con- 
sequently, the integrated spectra of m are inconsistent in 
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FIG. 7: Power-law exponents of Run III. See the caption of 
Fig. [6]for details. 



our simulations, and are in fact accumulation-dependent. 

In actual fact, the integrated spectra as functions the 
vertical wavenumbers are less steep than the GM spec- 
trum. Meanwhile, the cross-sectional spectra as func- 
tions the vertical wavenumbers in the inertial subrange 
are steeper than the GM spectrum. It comes from the 
fact that the exponents of the integrated spectra are a 
weighted mean of those of the cross-sectional spectra. 
In other words, /3int strongly depends on /3(1) and /3(2). 
Note that the horizontal wavenumbers k — corresponds 
to the inertial frequencies. 

The energy spectrum does not have double-power laws 
in the large wavenumbers in Run III. It can be explained 
by contamination of the inertial subrange in the small 
horizontal wavenumbers and the large density wavenum- 
bers, where the cross-sectional spectra do not show self- 
similarity in Fig. [5J It is curious that another double- 
power law appears in Run III in the small horizontal 
and density wavenumbers. Indeed, the spectrum at small 
wavenumbers can be approximated by 



E{k, m) oc k ^|m|" 



(30) 



As this power law does not appear in the inertial sub- 
range, it has limited relevance to the scope of this paper. 



FIG. 8: Power-law exponents of each cross-sectional spec- 
trum of the kinetic energy in Run II. left: A(|m|), which 
are the exponents of cross-sectional energy spectra Km{i^) in 
Fig- El right: B{ijj), which are the exponents of cross-sectional 
energy spectra A'ij(|m|) in Fig. [2] 



Similar to the cross-sectional spectra in (fc, m) space, 
the spectral exponents in (w, m) space could also be mea- 
sured. The power-law exponents of the integrated and 
cross-sectional spectra of the kinetic energy in Figs.[5]and 
[3] are shown in Figs. [5] and O The least-square fittings 
are made in cj € [10~^, 10~'^]rad/sec and |m| £ [8,32] 
as Km{^) oc tj-^dml) g^j_^(j K^{m) cx m-^d^D. We cannot 
find any characteristic exponents in the figures. There- 
fore it appears that the spectral energy density of these 
runs is not separable in the (w, m) space. This statement 
appears to be consistent with the fact that the cross- 
sectional spectra of the kinetic energy are neither self- 
similar nor separable. 

To summarize, Runs II and III exhibit the accumula- 
tion of energy at the near-inertial waves. The formation 
mechanism of the near-inertial accumulation is consis- 
tent with the oceanographic PSI arguments. Further- 
more, the cross-sectional spectra in Figs. |4] and [5] appear 
to be more separable than those in Figs. [5] and [31 It is 
supported by the fact that constant-exponent regions ap- 
pear in the spectra in (fc, m) space (Figs. [51 and ^ but 
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FIG. 9: Power-law exponents of each cross-sectional spec- 
trum of the kinetic energy in Run III. See the caption of 
Fig. [8] for details. 



not in in the spectra in {uu.m) space (Figs. [S] and 0- 
Therefore, the spectrum appear to be more separable in 
{k,m) space than in (a;,m) space in Run II and III. It 
also appears that the spectrum in the inertial subrange is 
sensitive to the details of the near-inertial accumulation. 
Most important, we reproduce the —2 exponent of the 
integrated spectra of the horizontal wavcnumbers. 



C. Run IV 

To further illustrate the influence of the accumulation 
around the horizontally longest waves, and to investi- 
gate the importance of the nonlocal interactions in the 
wavenumber space, we perform Run IV of TablelH There 
we choose the initial condition same as the Run I but with 
no energy in fc < 3 and |m| < 16. The energy spectrum 
after about 1.5 days is denoted by En\{k, \m\). The en- 
ergy spectrum developed from the GM spectrum in Run I 
that is shown in Fig. [If right) is denoted by Ecuik, |m|). 
Figure [10] shows the relative difference defined as 
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(31) 



FIG. 10: Relative difference between energy spectra after 
about 1.5 days developed from the GM spectrum with and 
without energy in small horizontal wavenumbers -Bd(fc, \m\). 



The wavenumbers different in the initial conditions ap- 
pear as a black rectangle in bottom left in Fig. [TOl If 
the nonlocal interactions with the accumulation of en- 
ergy around the horizontally longest waves were not dom- 
inant, the relative difference in the inertial subrange in 
Fig. [TOl would be small or slightly negative since the non- 
linear interactions try to compensate the defect in the 
small wavenumbers. Instead, E^i^k, \m\) has more than 
20% energy in 10 < \m\ < 40 and less energy transfer 
to the dissipation region in |m| > 80. This suggests that 
the small wavenumbers in the accumulation of energy 
transfer energy from small density wavenumbers to large 
density wavenumbers in the inertial subrange. This sce- 
nario of the energy transfer to the large density wavenum- 
bers is qualitatively consistent with the Induced Diffusion 
mechanism^i^^. 



D. Run V 

Run II and III have indicated that the external forcing 
at multiple frequencies does not produce smooth spec- 
tral energy density in the frequency space. We therefore 
are led to the question what would happen if the forced 
wavenumbers correspond to the single frequency. In Run 
V, the forced wavenumbers are modeled as M2 tides and 
have an almost single frequency of 3/. The single fre- 
quency forcing is significantly different from the band of 
forced frequencies in Run II and III. One of the reasons 
to choose 3/ is to further confirm that PSI mechanism is 
effective. Indeed, PSI is dominant in transferring energy 
to small k and large m wavenumbers when the forced 
frequencies are greater than 2/. Consequently, with 3/ 
forcing PSI will create a second peak at 1.5/. Although 
one could choose 2/ as a forced frequency, we have chosen 
3/ in order to be able to distinguish the direct excitations 
that appear at 1.5/ via PSI and the accumulation of en- 
ergy at the near-inertial frequencies around /. 

Fixing forced frequency to be 3/, we still have a free- 
dom of choosing forced wavenumbers. For simplicity, to 
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FIG. 11: Two-dimensional energy spectrum E{k, |m|) of Run 
V. The contours are plotted every powers of ten from lO"^'^ to 
1. The straight line from the origin shows the wavenumbers 
uj — 1.5/. The inset is the enlargement of the area near the 
origin. 



implement the forcing numerically, the external forcing 
is added to 24 small wavenumbers whose frequencies are 
close to 3/. Specifically, we keep the amplitudes of canon- 
ical variables, a(p), to be constant in time for the fol- 
lowing wavenumbers: {k,m) — (±1,0, ±2), (0,±1,±2), 
(±1,±1,±3), (±2,0, ±4) and (0,±2,±4). 

In Fig. [Tl] we show the two-dimensional energy spec- 
trum close to a statistically steady state obtained in this 
simulation. We observe that there is a significant energy 
accumulation around k 1, i.e. the horizontally longest 
waves. The wavenumbers that have 1 < fc < 2 except 
the forced wavenumbers, which are considered as the ac- 
cumulation of energy, account for approximately 60% of 
the total energy. Appearance of the accumulation in the 
near-inertial waves is consistent with oceanographic ob- 
servations. It is also consistent with Runs II and III. This 
consistency in observing the accumulation of energy in 
the near-inertial waves is one of the main results of the 
present paper. 

Aside from the accumulation of energy around the hor- 
izontally longest waves, other excited wavenumbers ap- 
pear around the line w = 1.5/. The frequency is half of 
the frequency of the forced wavenumbers. This excita- 
tion can be qualitatively be explained by the PSI mech- 
anism. Indeed, the excited wavenumbers whose frequen- 
cies are close to 1.5/ make resonant triads with the forced 
wavenumbers via PSI. Furthermore, another excitation 
can be seen around m ~ 5. The excitation can be inter- 
preted as resonant interactions with the accumulation of 
energy due to the ES mechanism. 

Figure [T^ shows integrated spectra, Eint{k) and 
£^int(|'7i|), and cross-sectional spectra, Em{k) and 
-Efe(|m|), obtained from the two-dimensional energy spec- 
trum shown in Fig. 1111 The integrated spectra can be 
interpreted roughly to be 




I 2 4 a 1C 32 M 1^ 



FIG. 12: Integrated and cross-sectional spectra of Run V. 
See the caption of Fig. |4]for details. 



and 

:Bi„t(|m|)cx |m|-i-33±0-3i. 

The exponents are obtained with the least-square method 
in the intervals of fc S [8,128] for £'int(fc) and that of 
|m| e [4,32] for i?int(|"i|). Note that the exponents 
are not too far from the large-wavenumber self-similar 
form of the GM spectrum. Indeed, the large wavenum- 
ber asymptotic form of the GM spectrum is given by 
Eq. 

We emphasize nevertheless that similarity of the in- 
tegrated numerical spectrum as a function of vertical 
wavenumbers with Eq. (|26|) may be coincidental to some 
degree. The exponent of the integrated spectrum of the 
GM spectrum in large vertical wavenumbers is not —I 
but —2. In addition, as apparent from Fig.[T^ the spec- 
trum in (fc, m) space is not separable. Moreover, as ex- 
plained above, the exponent of the vertical wavenum- 
bers is mostly determined by the accumulation of energy 
at small fc. Consequently it is rather sensitive to the 
specifics of numerical experiments. On the other hand, 
the fc exponent is relatively insensitive to the accumu- 
lation and is determined by the inertial subrange. The 
behavior was also observed in Run II and III. 
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FIG. 13: Integrated kinetic energy spectra, and cross- 
sectional energy spectra of Run V. The cross-sectional spec- 
tra that are functions of frequencies (left) are shown every 
eight curves for visibility. The cross-sectional spectra that 
are functions of vertical wavenumbers (right) are shown when 
LO = 2.5, 5, 10, 20, 40 x lO'^rad/sec. See the caption of Fig. [2] 
for details. 



We also measure the kinetic energy spectrum in (w, m) 
space. Figure [T51 shows the results for Run V. It appears 
that some energy exists in the large frequencies greater 
than the buoyancy frequency. It comes mainly from the 
non-pcriodicity of the time series, and is called the Gibbs 
phenomenon. 

Observe that the integrated spectrum of frequency is 
fitted by 



i?int(w) OC LO' 



-2.19±0.11 



over u! G [10~^, 10"'^]. The asymptotic behavior is 
consistent with the GM spectrum and oceanographic 
observations. This is the main result of this paper. 
Also note that the vertical spectrum could be fitted by 

|^|-1.12±0.35 oygj. g [4^32], 

Furthermore, one can interpret Fig. [13] as having three 
peaks in i?int(w). The accumulation of energy makes 
the peak around the inertial frequency / of this run, 
namely 5 x 10~^rad/sec. The second peak can be found 
around 7 x 10~^rad/sec which is the half of the frequency 



FIG. 14: Power-law exponents of Run V. See the caption of 
Fig. [S]for details. 



of the forced wavenumbers. Appearance of the exci- 
tation at half the forcing frequency is consistent with 
the PSI mechanism. The third peak is seen around 
1.4 X 10~''rad/sec which is the frequency of the forced 
wavenumbers. 

To illustrate the spectral details in the inertial sub- 
range, we show the power-law exponents of each cross- 
sectional spectrum in Fig. [TU a(|m|) and P{k). The 
exponents are obtained by fitting Em{k) (x in 
k e [16,128] and Ek{m) oc jmpC^) in |m| G [16,64]. The 
exponents of the integrated spectra, aint — —1.98 and 
/3i„t = —1.33, are also shown in Fig. 1141 The main point 
taken from the figure is that the integrated spectra have 
different exponents from the ones of the cross-sectional 
spectra in the inertial subrange. We observe that a varies 
between —1 and —5. Similarly, (3 varies between —1.5 and 
-3.5. 

Similar to Runs 11 and III, the exponent of the inte- 
grated spectrum of the horizontal wavenumbers is close to 
the GM spectrum. However, the cross-sectional spectra 
are steeper than the GM spectrum in the large horizon- 
tal and vertical wavenumbers considered as the inertial 
subrange. This can be explained by the fact that they 
are contaminated by direct excitations due to PSI and 
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FIG. 15: Power-law exponents of the kinetic energy in Run 
V. See the caption of Fig. [8] for details. 



ES mechanisms. Apparently, the cross-sectional expo- 
nents are not constant. Therefore, the energy spectrum 
of the numerical simulation cannot be accurately fitted 
by double-power functions. 

As the integrated and cross-sectional spectra of the 
GM spectrum have different exponents with respect to 
the vertical wavenumbers, the exponents of the inte- 
grated spectra are much different those of the cross- 
sectional spectra in the inertial subrange also in our 
simulations. This discrepancy is caused by the non- 
separability of the two-dimensional spectrum. The dis- 
crepancy in the vertical-wavenumber exponents is espe- 
cially clear by the near-inertial accumulation. The in- 
tegrated spectrum, Eintim), is established mainly by 
Ek=i,2{m) i.e. the accumulation of energy. Therefore, 
the exponent of the integrated spectrum is determined 
exclusively by the exponents of the cross-sectional spec- 
tra in the small horizontal wavenumbers. 

Figure [15] shows power-law exponents of each cross- 
sectional spectrum of the kinetic energy as functions of uj 
and \m\. The fits are made over uj G [10"'', 10~^]rad/sec 
or over \m\ £ [8,32]. The exponents of the integrated 
spectra are also shown in the figure for reference. It is in- 
structive to compare Figs.[Hland[T51 It appears that the 



(uj,m) spectrum is more separable than the (k,m) spec- 
trum. Indeed, the exponent a of horizontal wavenumbers 
of the {k,in) spectrum of Fig. [M] varies between —1 and 
—5, while the exponent A of frequencies of the (w, m) 
spectrum of Fig. 1151 varies between —1.5 and —4. Fur- 
thermore, it appears that A can loosely be represented 
as having the value of ^ ~ —2.7. Similarly, apart from 
the small- wavenumber and low-frequency region, the ex- 
ponent (3 of vertical wavenumbers of the (fc, m) spectrum 
varies between —2.5 and —3.5, and the exponent B of 
vertical wavenumbers of the (w, m) spectrum can loosely 
be interpreted as having value around —3.5. However, 
the separability ([T]) can not be interpreted as being fully 
satisfied. Our observation of the spectrum of Run V be- 
ing more separable in (w,rn) space than in (fc,m) space 
is consistent with observationally-based intuition. We 
again note that the exponents of the cross-sectional spec- 
tra are also far from those of the integrated spectra. This 
statement is true for both {k,m) and (w,m) spectra. 

As explained above, the near-inertial accumulation ap- 
pears only in the small horizontal wavenumbers and the 
small frequencies. Thus, the accumulation does not 
affect the integrated spectrum of k in the large hori- 
zontal wavenumbers when the integration of the two- 
dimensional {k,m) spectrum over m to obtain the inte- 
grated spectrum is made. Similarly, the integrated spec- 
trum of Lu in high frequencies are not affected by the 
near-inertial accumulation. Then, we can argue that the 
integrated spectra of the horizontal wavenumbers and the 
frequencies are determined by the wavenumbers in the 
inertial subrange. Consequently these integrated spectra 
appear to be characteristic in our numerical runs. There- 
fore, the power-law exponents of the integrated spectra 
roughly corresponding to those of the observations in- 
cluding the GM spectrum reflects the nonlinear interac- 
tions of the internal- wave system. 

To summarize. Run V is the largest and longest run 
that we have performed. The internal-wave field is 
forced by single frequency, small horizontal and verti- 
cal wavenumbers. Run V clearly demonstrates the en- 
ergy accumulation around the horizontally longest waves. 
This run also shows that PSI mechanism is effective at 
creating a second peak at half the forced frequency. Fur- 
thermore, this run reproduce the integrated spectra, k~'^ 
spectrum and w"^ spectrum. This is consistent with ob- 
servations. This run also demonstrates that the spectrum 
in {uj,m) space is more separable than in {k,m) space. 
It is also consistent with oceanic observations. 



V. DISCUSSION 

Stratified rotating turbulence is a complicated and fas- 
cinating subject. It has been a subject of intensive re- 
search in the last few decades. In particular various 
numerical simulations were performed. For example, 
Winters and D'Asaro^S, numerically modeled stratified 
Navier-Stokes turbulence. They consider a depth depen- 



16 



dent buoyancy frequency, and focus on a depth depen- 
dence of various physical quantities. Their simulation is 
done on a 32 X 256 x 129 grid. They implemented real- 
istic boundary conditions in the vertical directions and 
obtained the energy dissipation rates consistent with the- 
ories and observations. In contrast, our numerical simu- 
lations focus on spectral energy density. 

More recent simulations were also performed. Laval 
et al.— studied numerically stratified turbulence, and 
showed formation of pancake vortices and considerable 
amount of energy in the horizontally uniform vertical 
shear and vortical modes. Numerical simulations were 
also performed by Smith and Waleffeii for rotating strat- 
ified turbulence and by Smith and Lee^ for rotating tur- 
bulence. These simulations demonstrate the tendency 
of energy to accumulate at the horizontally large-scale 
flows. They presented that the accumulation is caused 
mainly by resonant wave-wave interactions. They also 
obtained one-dimensional spectra similar to our inte- 
grated spectra. 

In contrast to most of these simulations, our simula- 
tions are performed in periodic boxes in all three direc- 
tions. Furthermore, we completely exclude vortices and 
horizontally uniform vertical shears. In addition, we use 
larger numerical grids, which are required in order to an- 
alyze the inner structures of the energy spectra in the 
inertial subrange. In addition, larger numerical grids are 
necessary to make the nonlocal interactions more pro- 
nounced in the simulations. 

Certainly our reduced model constitutes a simplifica- 
tion of the ocean. One could argue that the horizon- 
tally uniform vertical shears excluded in our simulations 
might play an important role for the wavenumbers in 
the inertial subrange. Indeed, the horizontally uniform 
shear might strongly affect the spectra in the inertial sub- 
range as the accumulation of energy at the horizontally 
longest waves does. Furthermore, the largest-scale mo- 
tions also depend on other effects which are excluded 
from our simulation. In particular, we excluded variabil- 
ity of the inertial frequencies or (3 effect, seasonal vari- 
ability, and terrain properties. They may significantly 
affect the accumulation, providing mechanisms to dimin- 
ish or build the accumulation. Consequently the waves 
in the inertial subrange will be affected. 

As mentioned in Sec. |T1 anisotropic wave turbulence 
systems are often dominated by the nonlocal interactions. 
The oceanic internal waves are not an exception. There- 
fore it is important to study the effects of the nonlocal 
interactions in the wavenumber space. 

The invariance under Galilean transformation is bro- 
ken in the systems which have dominance of the nonlocal 
interactions. Indeed, Coriolis effect and presence of ver- 
tical stratification breaks the Galilean invariance of the 
internal waves. Consequently, one could argue that the 
dominance of the nonlocal interactions is caused by vi- 
olation of Galilean invariance of the systems. In fact, 
in the systems where Galilean invariance is broken, the 
large-scale flows not only advect the small-scale flows but 



also actively exchange energy with them almost like an 
inertia force. It is in contrast with the typical behavior 
of the large-scale motions in Galilean- invariant systems. 
In Galilean-invariant systems the large-scale flows advect 
small-scale flows, as in the sweeping in Navier-Stokes tur- 
bulence. 



VI. SUMMARY 

In this work, we performed numerical modeling of the 
reduced wave-only stratified rotating turbulence in order 
to investigate the wave-wave interactions. In particular, 
we made a series of 5 numerical runs with different initial 
conditions, grid sizes and forcing. 

Run I is the freely decaying run with the GM spectrum 
employed as the initial condition. It has demonstrated 
that the GM spectrum is not a steady-state spectrum for 
our numerical model. 

Runs II and III are forced-damped runs with different 
values of /. These runs consistently demonstrated the 
tendency of the energy to accumulate in the horizontally 
longest (i.e. near-inertial) waves. In these runs the level 
of the accumulation depends on the relation between the 
value of the inertial frequency / and forced frequencies. 
The accumulation of energy in the near-inertial frequen- 
cies is characteristic of the oceans. The runs also demon- 
strated largely non-separable spectra. It also appears 
that the energy spectrum in {k,m) space is more sepa- 
rable than in (u;,m) space. This is at odds with oceanic 
intuition. We also observed realistic k^'^ spectrum in 
Run II. Furthermore, the formation of the inertial peak 
is qualitatively consistent with a "named" nonlocal inter- 
action, which is the Parametric Subharmonic Instability. 
It also appears that the integrated vertical-wavenumber 
spectrum is highly sensitive to the accumulation at the 
horizontally longest waves. Finally, the integrated fre- 
quency spectra cannot be fitted by power laws due to 
the forcing at multiple frequencies. 

Run IV is the freely decaying run with the GM spec- 
trum without the small wavenumbers as initial condi- 
tions. This run demonstrates that the internal- wave sta- 
tistical properties are strongly affected by the nonlocal 
interactions with the small wavenumbers. In particular, 
removing the small wavenumbers from the initial con- 
dition produces a significantly different outcome. Thus, 
the energy transfer in the inertial subrange is qualita- 
tively explained by the Induced Diffusion. 

Run V also demonstrates that energy accumulates at 
the horizontally longest near-inertial waves. This is con- 
sistent with Runs II and III as well as the ocean. It also 
appears that the two-dimensional spectrum is more sepa- 
rable in {lu, m) space than in (fc, to) space. This is consis- 
tent with observationally-based intuition. Furthermore, 
the integrated frequency spectrum of Run V is given by 
Lu~'^, consistent with the ocean. Similarly, the integrated 
horizontal-wavenumber spectrum is also consistent 
with the observations. However, the integrated vertical- 
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wavenumber spectrum do not show the power-law ex- 
ponent consistent with the observations. This could 
be explained by the fact that the integrated vertical- 
wavenumber spectrum is determined by the accumula- 
tion at the near-inertial waves. The near-inertial waves 
may be much affected by processes not considered in this 
manuscript. 

In short, our numerical runs reproduce the following 
behavior: 

• Energy tends to accumulate at the horizontally 
longest waves i.e. near-inertial waves 

• To a lesser degree, some energy was observed in the 
small vertical wavenumbers. 

• Spectra in the inertial subrange are not completely 
separable 

• Realistic integrated horizontal k^^ spectrum is ob- 
served 

• Realistic integrated frequency spectrum is ob- 
served 

• The integrated vertical-wavenumber spectrum is 
determined exclusively by the near-inertial waves 

• The spectrum in the inertial subrange is largely de- 
termined by the interactions with the near-inertial 
waves. 



• The accumulation of energy of the near-inertial 
waves and the spectrum in the inertial subrange 
could qualitatively be described by "named" non- 
local interactions. 

The stratified rotating turbulence governed by the 
Navier-Stokes equation and its relation with the ocean 
will of course be further researched and debated. Here 
we address the reduced wave-wave numerical model of 
stratified rotating turbulence with the hope that it could 
shed light on processes in oceans and its spectral en- 
ergy density in particular. It appears that this reduced 
model reproduce key features of spectral energy density 
of oceanic internal waves. This supports observationally- 
based intuition that the spectral energy density of inter- 
nal waves is described predominantly by the wave-wave 
interactions. The future will certainly bring many further 
exciting developments, and the synthesis of theoretical, 
observational and numerical results yet to be obtained. 
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